Datos de estudio
En este ejemplo práctico utilizaremos datos de un estudio sobre la expresión génica en células infectadas con el virus de la gripe (Bajwa et al., 2016). En este estudio se infectó células dendríticas plasmocitoides humanas con el virus de la gripe y se comparó la expresión génica en estas células con la expresión génica en células no infectadas. El objetivo era ver cómo el virus de la gripe afectaba la función de estas células del sistema inmunitario.
Organización de las prácticas
En esta práctica vamos a usar el entorno de desarrollo integrado RStudio
(IDE) para R.
R es un lenguaje de programación para computación estadística y
gráficos.
El documento actual en el que estamos trabajando es un documento R Markdown. Son totalmente reproducibles y permiten combinar texto, imágenes y código —en lenguaje R.
Para renderizar un documento R Markdown a un archivo
HTML, solo necesitas hacer clic en el botón Knit que verás
en la barra de RStudio. Este archivo HTML puede compartirse como un
informe.
No necesitas renderizar todo el documento cada vez que quieras ver el
resultado de tu código en R. Puedes hacer clic en el botón
Run current chunk o usar el atajo de teclado
Ctrl+Alt+C, y el resultado del
código aparecerá debajo.
A lo largo del documento verás diferentes iconos, cuyo significado es:
: información adicional o útil
: un ejemplo
: un ejercicio práctico
: un espacio para responder al
ejercicio
: una pista para resolver un
ejercicio
Instalación de herramientas
Recomendamos encarecidamente usar un sistema operativo
Linux y acostumbrarse a trabajar con la
terminal.
En esta práctica, sin embargo, no la usaremos y casi no hay diferencias
entre usar RStudio en Windows o en Linux.
Instalación de R
El lenguaje de programación R puede descargarse desde este enlace, y está disponible para Windows, Linux y macOS.
Instalación de RStudio
RStudio se puede descargar desde este
enlace.
Se instala fácilmente en Windows, Linux y macOS.
Instalación de un paquete
Los paquetes de R son colecciones de funciones y/o datos desarrollados por la comunidad.
Para instalar un paquete, usamos la función
install.packages(), indicando entre comillas el nombre del
paquete que queremos instalar.
Instalación de ggplot2
ggplot2 es un paquete de visualización de datos para el
lenguaje de programación R (Wickham
2009).
Fue creado por Hadley Wickham, implementando la Grammar of
Graphics de Leland Wilkinson —un esquema general para la
visualización de datos que divide los gráficos en componentes semánticos
(Wilkinson
2010).
Para instalar el paquete ggplot2, utiliza lo
siguiente:
install.packages("ggplot2")
Cargar un paquete
Para cargar un paquete, usamos la función library(),
indicando el nombre de la librería que queremos cargar.
En el caso de ggplot2, usarás:
library("ggplot2")
¡Ya estás listo para empezar a usar R!
Introducción a R y RStudio
Directorio de trabajo
El directorio de trabajo es la carpeta de tu ordenador donde R busca por defecto los archivos que necesitas importar y donde guardará los resultados que exportes (tablas, gráficos, etc.).
Para saber cuál es el directorio de trabajo
(working directory) utilizamos la funcion
getwd():
## [1] "C:/Users/Marta/OneDrive - UAB/Bioestadística/Prácticas/P1"
Podemos definir una nueva ruta de trabajo con
setwd():
O de manera gráfica dentro de RStudio a Session > Set Working Directory > Choose Directory…).
Actividad 1.1
Decide una ruta de trabajo en tu ordenador con
setwd() o de manera gráfica. Comprueba que sea correcta con
getwd().
Por ejemplo, puedes escoger la carpeta donde te hayas descargado el material del Moodle.
En este directorio de trabajo que hayas escogido, realizaremos todos los análisis.
Leer el archivo con los datos de expresión
Lee en R el archivo “GSE68849_expression.csv”
que encontrarás en la carpeta data.
Puedes abrirlo primero con Excel o con otro programa de visualización de
datos (como un editor de texto) para comprobar si tiene nombres de
columnas y cómo están separadas.
Actividad 1.2
¿Cuál será el código en R para abrir el archivo?
Guarda el contenido del archivo en un objeto llamado
expressionData.
Puedes utilizar los comandos que vimos en
la diapositiva 12 de teoría del Tema 2: Introducción a R. Ten en cuenta
que quizás necesites cargar una librería, como
library(readr).
Actividad 1.3
Explora el archivo. Utiliza comandos como
head(expressionData) y View(expressionData)
para explorar el contenido de expressionData. ¿Qué
información contiene?
## # A tibble: 6 x 4
## subject treatment gene expression_level
## <chr> <chr> <chr> <dbl>
## 1 GSM1684095 control IFNA5 83.1
## 2 GSM1684096 influenza IFNA5 10096.
## 3 GSM1684097 control IFNA5 97.8
## 4 GSM1684098 influenza IFNA5 8181.
## 5 GSM1684099 control IFNA5 81.7
## 6 GSM1684100 influenza IFNA5 7055.
Tipos de variables
Actividad 1.4
Para saber el tipo de variable que contiene
expressionData, podemos utilizar la función
str(). Rellena la siguiente tabla indicando qué tipo de
variable contiene expressionData.
| Variable | Tipo de variable (Cuantitativa / Cualitativa) | Tipo específico (Discreta / Continua / Nominal / Ordinal) | Observaciones |
|---|---|---|---|
| Subject | Cualitativa | Nominal | Identificador único de cada sujeto |
| Treatment | Cualitativa | Nominal | Dos categorías: “control” e “influenza” |
| Gene | Cualitativa | Nominal | Nombre del gen analizado |
| Expression level | Cuantitativa | Continua | Niveles de expresión génica |
Estadística descriptiva
R contiene muchas funciones built-in para realizar
estadística descriptiva básica.
Actividad 1.5
Dada la columna con el nivel de expresión de los genes
expressionData$expression_level.
Calcula la media, la mediana y la
desviación estándar del vector
expressionData$expression_level.
## [1] 3947.509
## [1] 131.8218
## [1] 7210.434
Tidyverse
Actividad 1.6
¿Cuántos genes distintos hay en el set de datos?
Primero, podemos seleccionar la columna gene de nuestro
dataframe expressionData y luego usar
distinct() para quedarnos solo con los valores
únicos.
Después, aplicamos count() para contar cuántos elementos
hay en total.
Como queremos contar todas las filas resultantes, no hace falta
especificar ningún argumento dentro de count().
Recuerda usar %>% para
concatenar comandos.
## # A tibble: 1 x 1
## n
## <int>
## 1 10
Alternativamente, te mostramos una función que no se ha visto en la
teoría y que permite hacer lo mismo: n_distinct junto con
summarize:
## # A tibble: 1 x 1
## n
## <int>
## 1 10
Actividad 1.7
Una de las columnas es subject. Extrae todos los datos
correspondientes al sujeto GSM1684096 y guárdalos en un
nuevo objeto llamado expressionData_GSM1684096:
Utiliza la función
filter()
Actividad 1.8
Explora expressionData. Estamos
interesados en conocer el nivel de expresión del gen IFNA5 en
condición de infección para el sujeto GSM1684096. ¿Cómo lo
harías? Combina filter() y select():
expressionData_GSM1684096 %>% filter(gene=="IFNA5", treatment=="influenza") %>% select(expression_level)## # A tibble: 1 x 1
## expression_level
## <dbl>
## 1 10096.
Actividad 1.9
Calcula la media de expresión para cada tratamiento. Deberás combinar
group_by() y summarize(). La diapositiva 42
puede servir de guía:
## # A tibble: 2 x 2
## treatment avg
## <chr> <dbl>
## 1 control 156.
## 2 influenza 7739.
Actividad 1.10
Crea una nueva columna basada en el nivel de expresión de los genes: <10 → “bajo”,10-100 → “moderado”, 100-1000 → “alto”, >1000 → “muy alto”
En la diapositiva 32 viste cómo hacerlo con ifelse()
para dos grupos. Ahora lo haremos con case_when(), que
permite trabajar con múltiples condiciones. Un ejemplo que puede
servirte:
Actividad 1.11
En qué grupo de tratamiento estan los genes con categoría de
expresión muy alto? en los controles, en los de influenza,
o ambos?
## # A tibble: 1 x 1
## treatment
## <chr>
## 1 influenza
Visualización de datos
Actividad 1.12
Con esta actividad aprenderás a hacer paso a paso un heatmap.
Haremos un heatmap (con geom_tile) con el número de
genes (gene) por categoría de subject usando
los datos del archivo expressionData.
En la siguiente imagen puedes ver el gráfico final que haremos:
Paso 1. Visualizar los
datos
Primero, crearemos el objeto del heatmap con la función
ggplot, y generaremos un gráfico. Queremos un heatmap donde
se muestren los diferentes sujetos en el eje x
(subject) y los genes en el eje y
(gene), y el color de cada celda refleje cuánto se expresa
cada gen dentro de cada sujeto
(fill = log(expression_level)). Para ello, utilizaremos
geom_tile().
Modifica los ??? del siguiente código:
exp.heatmap <- ggplot(data = expressionData,
mapping = aes(x = subject,
y = gene,
fill = log(expression_level))) +
geom_tile()
exp.heatmapAunque es un buen comienzo, podemos mejorar varias cosas:
- Los ejes se podrían visualizar mejor, ya que no podemos leer los nombres de los sujetos.
- Estaría bien saber qué muestras corresponden a control y cuáles a infección.
- Usar otra escala de colores.
Paso 2. Mejorar el gráfico
Para mejorar los ejes, haremos lo siguiente:
- Utilizaremos una etiqueta más clara para el título del eje
x(por ejemplo, “Sujeto”). - Giraremos las etiquetas del eje
xcon un ángulo de 45 grados para poder leer los nombres de los sujetos. - Omitiremos el título del eje
y.
exp.heatmap <- exp.heatmap +
labs(x = "Subjecte", y = "") +
theme(axis.text.x = element_text(angle = 45, vjust = 0.5)) # Rotar les etiquetes de l'eix x
exp.heatmapA continuación, necesitamos separar las muestras de control de las
muestras infectadas con el virus de la gripe. Para hacerlo, utilizaremos
la función facet_grid() de ggplot, usando la
variable treatment para separar los dos grupos.
# Faceting
exp.heatmap <- exp.heatmap +
# facet_grid makes two panels, one for control, one for flu:
facet_grid(. ~ treatment, scales = "free_x", space = "free_x")
exp.heatmap Añade
scales = "free_x", space = "free_x" dentro de
facet_grid() para eliminar los espacios vacíos del gráfico.
¡Prueba a no añadirlos y observa qué pasa!
Y finalmente, podemos cambiar la escala de color; por ejemplo,
podríamos utilizar la escala de color
scale_fill_viridis():
Si has seguido todos los pasos, ¡ahora tu gráfico debería ser exactamente como el que se proporcionó al inicio!
¿Qué genes están sobreexpresados en el grupo infectado?
Respuesta:
IFNW1, IFNA5, IFNA2, IFNA16, IFNA13
Entrega
Entrega en el Campus Virtual el documento
renderizado (HTML) (botón knitr). Recuerda eliminar
eval = F de los chunks para poder visualizar el resultado
del código.
Referencias
Este ejemplo se ha creado
basándose en el tutorial
de Jeff Oliver para hacer el heatmap.